home *** CD-ROM | disk | FTP | other *** search
/ PD Collection CD 1 / PD Collection CD 1.iso / programer2 / icon / Source / Iconx / C / Rlargint < prev    next >
Encoding:
Text File  |  1990-07-20  |  45.3 KB  |  2,024 lines
  1. #include <math.h>
  2. #include "../h/config.h"
  3. #include "../h/rt.h"
  4. #include "rproto.h"
  5. #include <ctype.h>
  6.  
  7. #ifdef LargeInts
  8.  
  9. extern int over_flow;
  10.  
  11. /*
  12.  *  Conventions:
  13.  *
  14.  *  Bignums entering this module and leaving it are too large to
  15.  *  be represented with T_Integer.  So, externally, a given value
  16.  *  is always T_Integer or always T_Bignum.
  17.  *
  18.  *  Routines outside this module operate on bignums by calling
  19.  *  a routine like
  20.  *
  21.  *      bigadd(da, db, dx)
  22.  *
  23.  *  where da, db, and dx are pointers to tended descriptors.
  24.  *  For the common case where one argument is a T_Integer, these
  25.  *  call routines like
  26.  *
  27.  *      bigaddi(da, IntVal(*db), dx).
  28.  *
  29.  *  The bigxxxi routines can convert an integer to bignum form;
  30.  *  they use itobig.
  31.  *
  32.  *  The routines that actually do the work take (length, address)
  33.  *  pairs specifying unsigned base-B digit strings.  The sign handling
  34.  *  is done in the bigxxx routines.
  35.  */
  36.  
  37. /*
  38.  * Type for doing arithmetic on (2 * NB)-bit nonnegative numbers.
  39.  *  Normally unsigned but may be signed (with NB reduced appropriately)
  40.  *  if unsigned arithmetic is slow.
  41.  */
  42.  
  43. /* The bignum radix, B */
  44.  
  45. #define B            ((word)1 << NB)
  46.  
  47. /* Bignum digits in a word */
  48.  
  49. #define WORDLEN  (WordBits / NB + (WordBits % NB != 0))
  50.  
  51. /* size of a bignum block that will hold an integer */
  52.  
  53. #define INTBIGBLK  sizeof(struct b_bignum) + sizeof(DIGIT) * WORDLEN
  54.  
  55. /* lo(uword d) :            the low digit of a uword
  56.    hi(uword d) :            the rest, d is unsigned
  57.    signed_hi(uword d) :     the rest, d is signed
  58.    dbl(DIGIT a, DIGIT b) : the two-digit uword [a,b] */
  59.  
  60. #define lo(d)        ((d) & (B - 1))
  61. #define hi(d)        ((uword)(d) >> NB)
  62. #define dbl(a,b)     (((uword)(a) << NB) + (b))
  63.  
  64. #if ((-1) >> 1) < 0
  65. #define signed_hi(d) ((word)(d) >> NB)
  66. #else
  67. #define signbit      ((uword)1 << (WordBits - NB - 1))
  68. #define signed_hi(d) ((word)((((uword)(d) >> NB) ^ signbit) - signbit))
  69. #endif
  70.  
  71. /* BigNum(dptr dp) : the struct b_bignum pointed to by dp */
  72.  
  73. #define BigNum(dp)   ((struct b_bignum *)&BlkLoc(*dp)->bignumblk)
  74.  
  75. /* LEN(struct b_bignum *b) : number of significant digits */
  76.  
  77. #define LEN(b)       ((b)->lsd - (b)->msd + 1)
  78.  
  79. /* DIG(struct b_bignum *b, word i): pointer to ith most significant digit */
  80.  
  81. #define DIG(b,i)     (&(b)->digits[(b)->msd+(i)])
  82.  
  83. /* ceil, ln: ceil may be 1 too high in case ln is inaccurate */
  84.  
  85. #undef ceil
  86. #define ceil(x)      ((word)((x) + 1.01))
  87. #define ln(n)        (log((double)n))
  88.  
  89. /* determine the number of words needed for a bignum block with n digits */
  90.  
  91. #define BigNeed(n)   ( ((sizeof(struct b_bignum) + ((n) - 1) * sizeof(DIGIT)) \
  92.                + WordSize - 1) & -WordSize )
  93.  
  94. /* copied from rconv.c */
  95.  
  96. #if !EBCDIC
  97. #define tonum(c)     (isdigit(c) ? (c)-'0' : 10+(((c)|(040))-'a'))
  98. #else                    /* !EBCDIC */
  99. int tonum(c)
  100. int c;
  101. {
  102. #ifdef SASC
  103.    const static char *alphanum = "0123456789abcdefghijklmnopqrstuvwxyz" ;
  104.    char *where;
  105.  
  106.    where = memchr(alphanum, tolower(c), 36);
  107.    if (where == 0) return -1;
  108.    return where - alphanum;
  109. #else                    /* SASC */
  110.    if(isdigit(c)) return (c - '0');
  111.    if( (c | ' ') >= 'A' & (c | ' ') <= 'I') return( c - 'A' );
  112.    if( (c | ' ') >= 'J' & (c | ' ') <= 'R') return( (c - 'J') + 9 );
  113.    if( (c | ' ') >= 'S' & (c | ' ') <= 'Z') return( (c - 'S') + 18 );
  114.    return 0;
  115. #endif                    /* SASC */
  116.    }
  117. #endif                    /* !EBCDIC */
  118.  
  119. /* copied from oref.c */
  120.  
  121. #define RandVal (RanScale*(k_random=(RandA*(long)k_random+RandC)&0x7fffffffL))
  122.  
  123.  
  124. hidden int mkdesc    Params((struct b_bignum *x, dptr dx));
  125. hidden novalue itobig    Params((word i, struct b_bignum *x, dptr dx));
  126.  
  127. hidden novalue decout    Params((FILE *f, DIGIT *n, word l));
  128.  
  129. hidden int bigaddi    Params((dptr da, word i, dptr dx));
  130. hidden int bigsubi    Params((dptr da, word i, dptr dx));
  131. hidden int bigmuli    Params((dptr da, word i, dptr dx));
  132. hidden int bigdivi    Params((dptr da, word i, dptr dx));
  133. hidden int bigmodi    Params((dptr da, word i, dptr dx));
  134. hidden int bigpowi    Params((dptr da, word i, dptr dx));
  135. hidden int bigpowii    Params((word a, word i, dptr dx));
  136. hidden word bigcmpi    Params((dptr da, word i));
  137.  
  138. hidden DIGIT add1    Params((DIGIT *u, DIGIT *v, DIGIT *w, word n));
  139. hidden word sub1    Params((DIGIT *u, DIGIT *v, DIGIT *w, word n));
  140. hidden novalue mul1    Params((DIGIT *u, DIGIT *v, DIGIT *w, word n, word m));
  141. hidden novalue div1    
  142.         Params((DIGIT *a, DIGIT *b, DIGIT *q, DIGIT *r, word m, word n));
  143. hidden novalue compl1    Params((DIGIT *u, DIGIT *w, word n));
  144. hidden word cmp1    Params((DIGIT *u, DIGIT *v, word n));
  145. hidden DIGIT addi1    Params((DIGIT *u, word k, DIGIT *w, word n));
  146. hidden novalue subi1    Params((DIGIT *u, word k, DIGIT *w, word n));
  147. hidden DIGIT muli1    Params((DIGIT *u, word k, int c, DIGIT *w, word n));
  148. hidden DIGIT divi1    Params((DIGIT *u, word k, DIGIT *w, word n));
  149. hidden DIGIT shifti1    Params((DIGIT *u, word k, DIGIT c, DIGIT *w, word n));
  150. hidden word cmpi1    Params((DIGIT *u, word k, word n));
  151.  
  152. hidden novalue bdzero    Params((DIGIT *dest, word l));
  153. hidden novalue bdcopy    Params((DIGIT *src, DIGIT *dest, word l));
  154.  
  155. /*
  156.  * mkdesc -- put value into a descriptor
  157.  */
  158.  
  159. static int mkdesc(x, dx)
  160. struct b_bignum *x;
  161. dptr dx;
  162. {
  163.    word xlen, cmp;
  164.    static DIGIT maxword[WORDLEN] = { 1 << ((WordBits - 1) % NB) };
  165.  
  166.    /* suppress leading zero digits */
  167.  
  168.    while (x->msd != x->lsd && *DIG(x,0) == 0)
  169.       x->msd++;
  170.  
  171.    /* put it into a word if it fits, otherwise return the bignum */
  172.  
  173.    xlen = LEN(x);
  174.  
  175.    if (xlen < WORDLEN ||
  176.        (xlen == WORDLEN && ((cmp = cmp1(DIG(x,0), maxword, WORDLEN)) < 0 ||
  177.         (cmp == (word)0 && x->sign)))) {
  178.       word val = -(word)*DIG(x,0);
  179.       word i;
  180.  
  181.       for (i = x->msd; ++i <= x->lsd; )
  182.          val = (val << NB) - x->digits[i];
  183.       if (!x->sign)
  184.      val = -val;
  185.       dx->dword = D_Integer;
  186.       IntVal(*dx) = val;
  187.       }
  188.    else {
  189.       dx->dword = D_Bignum;
  190.       BlkLoc(*dx) = (union block *)x;
  191.       }
  192.    return Success;
  193. }
  194.  
  195. /*
  196.  *  i -> big
  197.  */
  198.  
  199. static novalue itobig(i, x, dx)
  200. word i;
  201. struct b_bignum *x;
  202. dptr dx;
  203. {
  204.    x->lsd = WORDLEN - 1;
  205.    x->msd = WORDLEN;
  206.    x->sign = 0;
  207.  
  208.    if (i == 0) {
  209.       x->msd--;
  210.       *DIG(x,0) = 0;
  211.       }
  212.    else if (i < 0) {
  213.       word d = lo(i);
  214.  
  215.       if (d != 0) {
  216.          d = B - d;
  217.          i += B;
  218.          }
  219.       i = - signed_hi(i);
  220.       x->msd--;
  221.       *DIG(x,0) = d;
  222.       x->sign = 1;
  223.       }
  224.             
  225.    while (i != 0) {
  226.       x->msd--;
  227.       *DIG(x,0) = lo(i);
  228.       i = hi(i);
  229.       }
  230.  
  231.    dx->dword = D_Bignum;
  232.    BlkLoc(*dx) = (union block *)x;
  233. }
  234.  
  235. /*
  236.  *  string -> bignum 
  237.  */
  238.  
  239. word bigradix(sign, r, s, dx)
  240. int sign;                      /* '-' or not */
  241. int r;                          /* radix 2 .. 36 */
  242. char *s;                        /* input string */
  243. dptr dx;                        /* output T_Integer or T_Bignum */
  244. {
  245.    struct b_bignum *b;
  246.    DIGIT *bd;
  247.    word len;
  248.    int c;
  249.  
  250.    if (r == 0)
  251.       return CvtFail;
  252.    len = ceil(strlen(s) * ln(r) / ln(B));
  253.    if (blkreq(BigNeed(len)) == Error)
  254.       return CvtFail;
  255.    b = alcbignum(len);
  256.    bd = DIG(b,0);
  257.  
  258.    bdzero(bd, len);
  259.  
  260.    if (r < 2 || r > 36)
  261.       return CvtFail;
  262.  
  263.    for (c = *s++; isalnum(c); c = *s++) {
  264.       c = tonum(c);
  265.       if (c >= r)
  266.      return CvtFail;
  267.       muli1(bd, (word)r, c, bd, len);
  268.       }
  269.  
  270.    while (isspace(c))
  271.       c = *s++;
  272.  
  273.    if (c != '\0')
  274.       return CvtFail;
  275.  
  276.    if (sign == '-')
  277.       b->sign = 1;
  278.  
  279.    /* put value into dx and return the type */
  280.  
  281.    (void)mkdesc(b, dx);
  282.    return Type(*dx);
  283. }
  284.  
  285. /*
  286.  *  bignum -> real
  287.  */
  288.  
  289. double bigtoreal(da)
  290. dptr da;
  291. {
  292.    word i;
  293.    double r = 0;
  294.    struct b_bignum *b = &BlkLoc(*da)->bignumblk;
  295.  
  296.    for (i = b->msd; i <= b->lsd; i++)
  297.       r = r * B + b->digits[i];
  298.  
  299.    return (b->sign ? -r : r);
  300. }
  301. /*
  302.  *  real -> bignum
  303.  */
  304.  
  305. int realtobig(da, dx)
  306. dptr da, dx;
  307. {
  308.    double x = BlkLoc(*da)->realblk.realval;
  309.    struct b_bignum *b;
  310.    word i, blen;
  311.    word d;
  312.  
  313.    if (x > 0.9999 * MinLong && x < 0.9999 * MaxLong) {
  314.       MakeInt((word)x, dx);
  315.       return Success;        /* got lucky; a simple integer suffices */
  316.       }
  317.  
  318.    x = x > 0 ? x : -x;
  319.    blen = ln(x) / ln(B) + 0.99;
  320.    for (i = 0; i < blen; i++)
  321.       x /= B;
  322.    if (x >= 1.0) {
  323.       x /= B;
  324.       blen += 1;
  325.       }
  326.  
  327.    if (blkreq(BigNeed(blen)) == Error)
  328.       return Error;
  329.    b = alcbignum(blen);
  330.    for (i = 0; i < blen; i++) {
  331.       d = (x *= B);
  332.       *DIG(b,i) = d;
  333.       x -= d;
  334.       }
  335.      
  336.    b->sign = x < 0;
  337.    return mkdesc(b, dx);
  338. }
  339.  
  340. /*
  341.  *  bignum -> string
  342.  */
  343.  
  344. int bigtos(da, dx)
  345. dptr da, dx;
  346. {
  347.    struct b_bignum *a, *temp;
  348.    word alen = LEN(BigNum(da));
  349.    word slen = ceil(alen * ln(B) / ln(10));
  350.    char *p, *q;
  351.  
  352.    if (strreq(slen) == Error || blkreq(BigNeed(alen)) == Error) 
  353.       return CvtFail;
  354.    a = BigNum(da);
  355.    temp = alcbignum(alen);
  356.    if (a->sign)
  357.       slen++;
  358.    q = alcstr("",slen);
  359.    bdcopy(DIG(a,0), DIG(temp,0), alen);
  360.    p = q += slen;
  361.    while (cmpi1(DIG(temp,0), (word)0, alen))
  362.       *--p = '0' + divi1(DIG(temp,0), (word)10, DIG(temp,0), alen);
  363.    if (a->sign)
  364.       *--p = '-';
  365.    StrLen(*dx) = q - p;
  366.    StrLoc(*dx) = p;
  367.    return Cvt;
  368. }
  369.  
  370. /*
  371.  *  bignum -> file 
  372.  */
  373.  
  374. novalue bigprint(f, da)
  375. FILE *f;
  376. dptr da;
  377. {
  378.    struct b_bignum *a, *temp;
  379.    word alen = LEN(BigNum(da));
  380.    word slen, dlen;
  381.  
  382.    slen = (BlkLoc(*da)->bignumblk.lsd - BlkLoc(*da)->bignumblk.msd + 1);
  383.    dlen = slen * NB * 0.3010299956639812;    /* 1 / log2(10) */
  384.    if (dlen > MaxDigits) {
  385.       fprintf(f, "integer(~%ld)",dlen - 2);    /* center estimate */
  386.       return;
  387.       }
  388.  
  389.    if (blkreq(BigNeed(alen)) == Error) {
  390.       fatalerr(0, NULL);        /* not worth passing this one back */
  391.       }
  392.    temp = alcbignum(alen);
  393.    a = BigNum(da);
  394.    bdcopy(DIG(a,0), DIG(temp,0), alen);
  395.    if (a->sign)
  396.       putc('-', f);
  397.    decout(f, DIG(temp,0), alen);
  398. }
  399.  
  400. static novalue decout(f, n, l)
  401. FILE *f;
  402. DIGIT *n;
  403. word l;
  404. {
  405.    word i = divi1(n, (word)10, n, l);
  406.  
  407.    if (cmpi1(n, (word)0, l))
  408.       decout(f, n, l);
  409.    putc('0' + i, f);
  410. }
  411.  
  412. /*
  413.  *  da -> dx
  414.  */
  415.  
  416. int cpbignum(da, dx)
  417. dptr da, dx;
  418. {
  419.    struct b_bignum *a;
  420.    word alen = LEN(BigNum(da));
  421.    struct b_bignum *x;
  422.  
  423.    if (blkreq(BigNeed(alen)) == Error)
  424.       return Error;
  425.    x = alcbignum(alen);
  426.    a = BigNum(da);
  427.    bdcopy(DIG(a,0), DIG(x,0), alen);
  428.    x->sign = a->sign;
  429.    return mkdesc(x, dx);
  430. }
  431.  
  432. /*
  433.  *  da + db -> dx
  434.  */
  435.  
  436. int bigadd(da, db, dx)
  437. dptr da, db;
  438. dptr dx;
  439. {
  440.    struct b_bignum *x, *a, *b;
  441.    word alen, blen;
  442.    word c;
  443.  
  444.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  445.       alen = LEN(BigNum(da));
  446.       blen = LEN(BigNum(db));
  447.       if (blkreq(BigNeed(alen > blen ? alen + 1 : blen + 1)) == Error)
  448.      return Error;
  449.       a = BigNum(da);
  450.       b = BigNum(db);
  451.       if (a->sign == b->sign) {
  452.          if (alen > blen) {
  453.             x = alcbignum(alen + 1);
  454.             c = add1(DIG(a,alen-blen), DIG(b,0), DIG(x,alen-blen+1), blen);
  455.             *DIG(x,0) = addi1(DIG(a,0), c, DIG(x,1), alen-blen);
  456.             }
  457.          else if (alen == blen) {
  458.             x = alcbignum(alen + 1);
  459.             *DIG(x,0) = add1(DIG(a,0), DIG(b,0), DIG(x,1), alen);
  460.             }
  461.          else {
  462.             x = alcbignum(blen + 1);
  463.             c = add1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen+1), alen);
  464.             *DIG(x,0) = addi1(DIG(b,0), c, DIG(x,1), blen-alen);
  465.             }
  466.          x->sign = a->sign;
  467.          }
  468.       else {
  469.          if (alen > blen) {
  470.             x = alcbignum(alen);
  471.             c = sub1(DIG(a,alen-blen), DIG(b,0), DIG(x,alen-blen), blen);
  472.             subi1(DIG(a,0), -c, DIG(x,0), alen-blen);
  473.             x->sign = a->sign;
  474.             }
  475.          else if (alen == blen) {
  476.             x = alcbignum(alen);
  477.             if (cmp1(DIG(a,0), DIG(b,0), alen) > 0) {
  478.                (void)sub1(DIG(a,0), DIG(b,0), DIG(x,0), alen);
  479.                x->sign = a->sign;
  480.                }
  481.             else {
  482.                (void)sub1(DIG(b,0), DIG(a,0), DIG(x,0), alen);
  483.                x->sign = b->sign;
  484.                }
  485.             }
  486.          else {
  487.             x = alcbignum(blen);
  488.             c = sub1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen), alen);
  489.             subi1(DIG(b,0), -c, DIG(x,0), blen-alen);
  490.             x->sign = b->sign;
  491.             }
  492.          }
  493.       return mkdesc(x, dx);
  494.       }
  495.    else if (Type(*da) == T_Bignum)    /* bignum + integer */
  496.       return bigaddi(da, IntVal(*db), dx);
  497.    else if (Type(*db) == T_Bignum)    /* integer + bignum */
  498.       return bigaddi(db, IntVal(*da), dx);
  499.    else {                             /* integer + integer */
  500.       struct descrip td;
  501.       char tdigits[INTBIGBLK];
  502.  
  503.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  504.       return bigaddi(&td, IntVal(*db), dx);
  505.       }
  506. }
  507.  
  508. /*
  509.  *  da - db -> dx
  510.  */ 
  511.  
  512. int bigsub(da, db, dx)
  513. dptr da, db, dx;
  514. {
  515.    struct descrip td;
  516.    char tdigits[INTBIGBLK];
  517.    struct b_bignum *a, *b, *x;
  518.    word alen, blen;
  519.    word c;
  520.  
  521.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  522.       alen = LEN(BigNum(da));
  523.       blen = LEN(BigNum(db));
  524.       if (blkreq(BigNeed(alen > blen ? alen + 1 : blen + 1)) == Error)
  525.      return Error;
  526.       a = BigNum(da);
  527.       b = BigNum(db);
  528.       if (a->sign != b->sign) {
  529.          if (alen > blen) {
  530.             x = alcbignum(alen + 1);
  531.             c = add1(DIG(a,alen-blen), DIG(b,0), DIG(x,alen-blen+1), blen);
  532.             *DIG(x,0) = addi1(DIG(a,0), c, DIG(x,1), alen-blen);
  533.             }
  534.          else if (alen == blen) {
  535.             x = alcbignum(alen + 1);
  536.             *DIG(x,0) = add1(DIG(a,0), DIG(b,0), DIG(x,1), alen);
  537.             }
  538.          else {
  539.             x = alcbignum(blen + 1);
  540.             c = add1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen+1), alen);
  541.             *DIG(x,0) = addi1(DIG(b,0), c, DIG(x,1), blen-alen);
  542.             }
  543.          x->sign = a->sign;
  544.          }
  545.       else {
  546.          if (alen > blen) {
  547.             x = alcbignum(alen);
  548.             c = sub1(DIG(a,alen-blen), DIG(b,0), DIG(x,alen-blen), blen);
  549.             subi1(DIG(a,0), -c, DIG(x,0), alen-blen);
  550.             x->sign = a->sign;
  551.             }
  552.          else if (alen == blen) {
  553.             x = alcbignum(alen);
  554.             if (cmp1(DIG(a,0), DIG(b,0), alen) > 0) {
  555.                (void)sub1(DIG(a,0), DIG(b,0), DIG(x,0), alen);
  556.                x->sign = a->sign;
  557.                }
  558.             else {
  559.                (void)sub1(DIG(b,0), DIG(a,0), DIG(x,0), alen);
  560.                x->sign = 1 ^ b->sign;
  561.                }
  562.             }
  563.          else {
  564.             x = alcbignum(blen);
  565.             c = sub1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen), alen);
  566.             subi1(DIG(b,0), -c, DIG(x,0), blen-alen);
  567.             x->sign = 1 ^ b->sign;
  568.             }
  569.          }
  570.       return mkdesc(x, dx);
  571.       }
  572.    else if (Type(*da) == T_Bignum)     /* bignum - integer */
  573.       return bigsubi(da, IntVal(*db), dx);
  574.    else if (Type(*db) == T_Bignum) {   /* integer - bignum */
  575.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  576.       alen = LEN(BigNum(&td));
  577.       blen = LEN(BigNum(db));
  578.       if (blkreq(BigNeed(alen > blen ? alen + 1 : blen + 1)) == Error)
  579.      return Error;
  580.       a = BigNum(&td);
  581.       b = BigNum(db);
  582.       if (a->sign != b->sign) {
  583.          if (alen == blen) {
  584.             x = alcbignum(alen + 1);
  585.             *DIG(x,0) = add1(DIG(a,0), DIG(b,0), DIG(x,1), alen);
  586.             }
  587.          else {
  588.             x = alcbignum(blen + 1);
  589.             c = add1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen+1), alen);
  590.             *DIG(x,0) = addi1(DIG(b,0), c, DIG(x,1), blen-alen);
  591.             }
  592.          x->sign = a->sign;
  593.          }
  594.       else {
  595.          if (alen == blen) {
  596.             x = alcbignum(alen);
  597.             if (cmp1(DIG(a,0), DIG(b,0), alen) > 0) {
  598.                (void)sub1(DIG(a,0), DIG(b,0), DIG(x,0), alen);
  599.                x->sign = a->sign;
  600.                }
  601.             else {
  602.                (void)sub1(DIG(b,0), DIG(a,0), DIG(x,0), alen);
  603.                x->sign = 1 ^ b->sign;
  604.                }
  605.             }
  606.          else {
  607.             x = alcbignum(blen);
  608.             c = sub1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen), alen);
  609.             subi1(DIG(b,0), -c, DIG(x,0), blen-alen);
  610.             x->sign = 1 ^ b->sign;
  611.             }
  612.          }
  613.       return mkdesc(x, dx);
  614.       }
  615.    else {                              /* integer - integer */
  616.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  617.       return bigsubi(&td, IntVal(*db), dx);
  618.       }
  619.       
  620. }
  621.  
  622. /*
  623.  *  da * db -> dx
  624.  */
  625.  
  626. int bigmul(da, db, dx)
  627. dptr da, db, dx;
  628. {
  629.    struct b_bignum *a, *b, *x;
  630.    word alen, blen;
  631.  
  632.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  633.       alen = LEN(BigNum(da));
  634.       blen = LEN(BigNum(db));
  635.       if (blkreq(BigNeed(alen + blen)) == Error)
  636.      return Error;
  637.       a = BigNum(da);
  638.       b = BigNum(db);
  639.       x = alcbignum(alen + blen);
  640.       mul1(DIG(a,0), DIG(b,0), DIG(x,0), alen, blen);
  641.       x->sign = a->sign ^ b->sign;
  642.       return mkdesc(x, dx);
  643.       }
  644.    else if (Type(*da) == T_Bignum)    /* bignum * integer */
  645.       return bigmuli(da, IntVal(*db), dx);
  646.    else if (Type(*db) == T_Bignum)    /* integer * bignum */
  647.       return bigmuli(db, IntVal(*da), dx);
  648.    else {                             /* integer * integer */
  649.       struct descrip td;
  650.       char tdigits[INTBIGBLK];
  651.  
  652.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  653.       return bigmuli(&td, IntVal(*db), dx);
  654.       }
  655. }
  656.  
  657. /*
  658.  *  da / db -> dx
  659.  */
  660.  
  661. int bigdiv(da, db, dx)
  662. dptr da, db, dx;
  663. {
  664.    struct b_bignum *a, *b, *x;
  665.    word alen, blen;
  666.  
  667.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  668.       alen = LEN(BigNum(da));
  669.       blen = LEN(BigNum(db));
  670.       if (alen < blen) {
  671.          MakeInt(0, dx);
  672.          return Success;
  673.          }
  674.       if (blkreq(BigNeed(alen-blen+1)+BigNeed(alen+1)+BigNeed(blen)) == Error)
  675.      return Error;
  676.       a = BigNum(da);
  677.       b = BigNum(db);
  678.       x = alcbignum(alen - blen + 1);
  679.       if (blen == 1)
  680.          divi1(DIG(a,0), (word)*DIG(b,0), DIG(x,0), alen);
  681.       else
  682.          div1(DIG(a,0), DIG(b,0), DIG(x,0), NULL, alen-blen, blen);
  683.       x->sign = a->sign ^ b->sign;
  684.       return mkdesc(x, dx);
  685.       }
  686.    else if (Type(*da) == T_Bignum)     /* bignum / integer */
  687.       return bigdivi(da, IntVal(*db), dx);
  688.    else if (Type(*db) == T_Bignum) {   /* integer / bignum */
  689.       MakeInt(0, dx);
  690.       return Success;
  691.       }
  692.    /* not called for integer / integer */
  693. }
  694.  
  695. /*
  696.  *  da % db -> dx
  697.  */
  698.  
  699. int bigmod(da, db, dx)
  700. dptr da, db, dx;
  701. {
  702.    struct b_bignum *a, *b, *x, *temp;
  703.    word alen, blen;
  704.  
  705.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  706.       alen = LEN(BigNum(da));
  707.       blen = LEN(BigNum(db));
  708.       if (alen < blen) {
  709.          cpbignum(da, dx);
  710.          return Success;
  711.          }
  712.       if (blkreq(BigNeed(blen)+BigNeed(alen+1)+BigNeed(blen)) == Error)
  713.      return Error;
  714.       a = BigNum(da);
  715.       b = BigNum(db);
  716.       x = alcbignum(blen);
  717.       if (blen == 1) {
  718.      temp = alcbignum(alen);
  719.          *DIG(x,0) = divi1(DIG(a,0), (word)*DIG(b,0), DIG(temp,0), alen);
  720.          }
  721.       else
  722.          div1(DIG(a,0), DIG(b,0), NULL, DIG(x,0), alen-blen, blen);
  723.       x->sign = a->sign;
  724.       return mkdesc(x, dx);
  725.       }
  726.    else if (Type(*da) == T_Bignum)     /* bignum % integer */
  727.       return bigmodi(da, IntVal(*db), dx);
  728.    else if (Type(*db) == T_Bignum) {   /* integer % bignum */
  729.       struct descrip td;
  730.       char tdigits[INTBIGBLK];
  731.  
  732.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  733.       cpbignum(&td, dx);
  734.       return Success;
  735.       }
  736.    /* not called for integer % integer */
  737. }
  738.  
  739. /*
  740.  *  -i -> dx
  741.  */
  742.  
  743. int bigneg(da, dx)
  744. dptr da, dx;
  745. {
  746.    struct descrip td;
  747.    char tdigits[INTBIGBLK];
  748.  
  749.    itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  750.    BigNum(&td)->sign ^= 1;
  751.    return cpbignum(&td, dx);
  752. }
  753.  
  754. /*
  755.  *  da ^ db -> dx
  756.  */
  757.  
  758. int bigpow(da, db, dx)
  759. dptr da, db, dx;
  760. {
  761.    word n;
  762.  
  763.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  764.       if (BigNum(db)->sign) {
  765.          MakeInt(0, dx);
  766.          }
  767.       else {
  768.          n = LEN(BigNum(db)) * NB;
  769.          /* scan bits left to right.  skip leading 1. */
  770.          while (--n >= 0)
  771.             if ((*DIG(BigNum(db), n / NB) & (1 << (n % NB))))
  772.                break;
  773.          /* then, for each zero, square the partial result;
  774.           *  for each one, square it and multiply it by a */
  775.          *dx = *da;
  776.          while (--n >= 0) {
  777.             if (bigmul(dx, dx, dx) == Error)
  778.            return Error;
  779.             if ((*DIG(BigNum(db), n / NB) & (1 << (n % NB))))
  780.                if (bigmul(dx, da, dx) == Error)
  781.           return Error;
  782.         }
  783.          }
  784.       return Success;
  785.       }
  786.    else if (Type(*da) == T_Bignum)    /* bignum ^ integer */
  787.       return bigpowi(da, IntVal(*db), dx);
  788.    else if (Type(*db) == T_Bignum)    /* integer ^ bignum */
  789.       return bigpowii(IntVal(*da), (word)bigtoreal(db), dx);
  790.    else                               /* integer ^ integer */
  791.       return bigpowii(IntVal(*da), IntVal(*db), dx);
  792. }
  793.  
  794. /*
  795.  *  iand(da, db) -> dx
  796.  */
  797.  
  798. int bigand(da, db, dx)
  799. dptr da, db, dx;
  800. {
  801.    struct b_bignum *a, *b, *x;
  802.    word alen, blen, xlen;
  803.    word i;
  804.    struct b_bignum *tad, *tbd;
  805.    DIGIT *ad, *bd;
  806.    struct descrip td;
  807.    char tdigits[INTBIGBLK];
  808.  
  809.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  810.       alen = LEN(BigNum(da));
  811.       blen = LEN(BigNum(db));
  812.       xlen = alen > blen ? alen : blen;
  813.       if (blkreq(3 * BigNeed(xlen)) == Error)
  814.          return Error;
  815.       a = BigNum(da);
  816.       b = BigNum(db);
  817.       x = alcbignum(xlen);
  818.  
  819.       if (alen == xlen && !a->sign)
  820.          ad = DIG(a,0);
  821.       else {
  822.          tad = alcbignum(xlen);
  823.          ad = DIG(tad,0);
  824.          bdzero(ad, xlen - alen);
  825.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  826.          if (a->sign)
  827.         compl1(ad, ad, xlen);
  828.          }
  829.  
  830.       if (blen == xlen && !b->sign)
  831.          bd = DIG(b,0);
  832.       else {
  833.          tbd = alcbignum(xlen);
  834.          bd = DIG(tbd,0);
  835.          bdzero(bd, xlen - blen);
  836.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  837.          if (b->sign)
  838.         compl1(bd, bd, xlen);
  839.          }
  840.         
  841.       for (i = 0; i < xlen; i++)
  842.          *DIG(x,i) = ad[i] & bd[i];
  843.  
  844.       if (a->sign & b->sign) {
  845.          x->sign = 1;
  846.          compl1(DIG(x,0), DIG(x,0), xlen);
  847.          }
  848.       }
  849.    else if (Type(*da) == T_Bignum) {   /* iand(bignum,integer) */
  850.       itobig(IntVal(*db), (struct b_bignum *)tdigits, &td);
  851.       alen = LEN(BigNum(da));
  852.       blen = LEN(BigNum(&td));
  853.       xlen = alen > blen ? alen : blen;
  854.       if (blkreq(3 * BigNeed(alen)) == Error)
  855.          return Error;
  856.       a = BigNum(da);
  857.       b = BigNum(&td);
  858.       x = alcbignum(alen);
  859.  
  860.       if (alen == xlen && !a->sign)
  861.          ad = DIG(a,0);
  862.       else {
  863.          tad = alcbignum(xlen);
  864.          ad = DIG(tad,0);
  865.          bdzero(ad, xlen - alen);
  866.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  867.          if (a->sign)
  868.         compl1(ad, ad, xlen);
  869.          }
  870.  
  871.       if (blen == xlen && !b->sign)
  872.          bd = DIG(b,0);
  873.       else {
  874.          tbd = alcbignum(xlen);
  875.          bd = DIG(tbd,0);
  876.          bdzero(bd, xlen - blen);
  877.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  878.          if (b->sign)
  879.         compl1(bd, bd, xlen);
  880.          }
  881.         
  882.       for (i = 0; i < xlen; i++)
  883.          *DIG(x,i) = ad[i] & bd[i];
  884.  
  885.       if (a->sign & b->sign) {
  886.          x->sign = 1;
  887.          compl1(DIG(x,0), DIG(x,0), xlen);
  888.          }
  889.       }
  890.    else if (Type(*db) == T_Bignum) {   /* iand(integer,bignum) */
  891.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  892.       alen = LEN(BigNum(&td));
  893.       blen = LEN(BigNum(db));
  894.       xlen = alen > blen ? alen : blen;
  895.       if (blkreq(3 * BigNeed(blen)) == Error)
  896.          return Error;
  897.       a = BigNum(&td);
  898.       b = BigNum(db);
  899.       x = alcbignum(blen);
  900.  
  901.       if (alen == xlen && !a->sign)
  902.          ad = DIG(a,0);
  903.       else {
  904.          tad = alcbignum(xlen);
  905.          ad = DIG(tad,0);
  906.          bdzero(ad, xlen - alen);
  907.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  908.          if (a->sign)
  909.         compl1(ad, ad, xlen);
  910.          }
  911.  
  912.       if (blen == xlen && !b->sign)
  913.          bd = DIG(b,0);
  914.       else {
  915.          tbd = alcbignum(xlen);
  916.          bd = DIG(tbd,0);
  917.          bdzero(bd, xlen - blen);
  918.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  919.          if (b->sign)
  920.         compl1(bd, bd, xlen);
  921.          }
  922.         
  923.       for (i = 0; i < xlen; i++)
  924.          *DIG(x,i) = ad[i] & bd[i];
  925.  
  926.       if (a->sign & b->sign) {
  927.          x->sign = 1;
  928.          compl1(DIG(x,0), DIG(x,0), xlen);
  929.          }
  930.       }
  931.    /* not called for iand(integer,integer) */
  932.  
  933.    return mkdesc(x, dx);
  934. }
  935.  
  936. /*
  937.  *  ior(da, db) -> dx
  938.  */
  939.  
  940. int bigor(da, db, dx)
  941. dptr da, db, dx;
  942. {
  943.    struct b_bignum *a, *b, *x;
  944.    word alen, blen, xlen;
  945.    word i;
  946.    struct b_bignum *tad, *tbd;
  947.    DIGIT *ad, *bd;
  948.    struct descrip td;
  949.    char tdigits[INTBIGBLK];
  950.  
  951.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  952.       alen = LEN(BigNum(da));
  953.       blen = LEN(BigNum(db));
  954.       xlen = alen > blen ? alen : blen;
  955.       if (blkreq(3 * BigNeed(xlen)) == Error)
  956.          return Error;
  957.       a = BigNum(da);
  958.       b = BigNum(db);
  959.       x = alcbignum(xlen);
  960.  
  961.       if (alen == xlen && !a->sign)
  962.          ad = DIG(a,0);
  963.       else {
  964.          tad = alcbignum(xlen);
  965.          ad = DIG(tad,0);
  966.          bdzero(ad, xlen - alen);
  967.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  968.          if (a->sign)
  969.         compl1(ad, ad, xlen);
  970.          }
  971.  
  972.       if (blen == xlen && !b->sign)
  973.          bd = DIG(b,0);
  974.       else {
  975.          tbd = alcbignum(xlen);
  976.          bd = DIG(tbd,0);
  977.          bdzero(bd, xlen - blen);
  978.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  979.          if (b->sign)
  980.         compl1(bd, bd, xlen);
  981.          }
  982.         
  983.       for (i = 0; i < xlen; i++)
  984.          *DIG(x,i) = ad[i] | bd[i];
  985.  
  986.       if (a->sign | b->sign) {
  987.          x->sign = 1;
  988.          compl1(DIG(x,0), DIG(x,0), xlen);
  989.          }
  990.       }
  991.    else if (Type(*da) == T_Bignum) {   /* ior(bignum,integer) */
  992.       itobig(IntVal(*db), (struct b_bignum *)tdigits, &td);
  993.       alen = LEN(BigNum(da));
  994.       blen = LEN(BigNum(&td));
  995.       xlen = alen > blen ? alen : blen;
  996.       if (blkreq(3 * BigNeed(alen)) == Error)
  997.          return Error;
  998.       a = BigNum(da);
  999.       b = BigNum(&td);
  1000.       x = alcbignum(alen);
  1001.  
  1002.       if (alen == xlen && !a->sign)
  1003.          ad = DIG(a,0);
  1004.       else {
  1005.          tad = alcbignum(xlen);
  1006.          ad = DIG(tad,0);
  1007.          bdzero(ad, xlen - alen);
  1008.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  1009.          if (a->sign)
  1010.         compl1(ad, ad, xlen);
  1011.          }
  1012.  
  1013.       if (blen == xlen && !b->sign)
  1014.          bd = DIG(b,0);
  1015.       else {
  1016.          tbd = alcbignum(xlen);
  1017.          bd = DIG(tbd,0);
  1018.          bdzero(bd, xlen - blen);
  1019.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  1020.          if (b->sign)
  1021.         compl1(bd, bd, xlen);
  1022.          }
  1023.         
  1024.       for (i = 0; i < xlen; i++)
  1025.          *DIG(x,i) = ad[i] | bd[i];
  1026.  
  1027.       if (a->sign | b->sign) {
  1028.          x->sign = 1;
  1029.          compl1(DIG(x,0), DIG(x,0), xlen);
  1030.          }
  1031.       }
  1032.    else if (Type(*db) == T_Bignum) {   /* ior(integer,bignym) */
  1033.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  1034.       alen = LEN(BigNum(&td));
  1035.       blen = LEN(BigNum(db));
  1036.       xlen = alen > blen ? alen : blen;
  1037.       if (blkreq(3 * BigNeed(blen)) == Error)
  1038.          return Error;
  1039.       a = BigNum(&td);
  1040.       b = BigNum(db);
  1041.       x = alcbignum(blen);
  1042.  
  1043.       if (alen == xlen && !a->sign)
  1044.          ad = DIG(a,0);
  1045.       else {
  1046.          tad = alcbignum(xlen);
  1047.          ad = DIG(tad,0);
  1048.          bdzero(ad, xlen - alen);
  1049.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  1050.          if (a->sign)
  1051.         compl1(ad, ad, xlen);
  1052.          }
  1053.  
  1054.       if (blen == xlen && !b->sign)
  1055.          bd = DIG(b,0);
  1056.       else {
  1057.          tbd = alcbignum(xlen);
  1058.          bd = DIG(tbd,0);
  1059.          bdzero(bd, xlen - blen);
  1060.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  1061.          if (b->sign)
  1062.         compl1(bd, bd, xlen);
  1063.          }
  1064.         
  1065.       for (i = 0; i < xlen; i++)
  1066.          *DIG(x,i) = ad[i] | bd[i];
  1067.  
  1068.       if (a->sign | b->sign) {
  1069.          x->sign = 1;
  1070.          compl1(DIG(x,0), DIG(x,0), xlen);
  1071.          }
  1072.       }
  1073.    /* not called for ior(integer,integer) */
  1074.  
  1075.    return mkdesc(x, dx);
  1076. }
  1077.  
  1078. /*
  1079.  *  xor(da, db) -> dx
  1080.  */
  1081.  
  1082. int bigxor(da, db, dx)
  1083. dptr da, db, dx;
  1084. {
  1085.    struct b_bignum *a, *b, *x;
  1086.    word alen, blen, xlen;
  1087.    word i;
  1088.    struct b_bignum *tad, *tbd;
  1089.    DIGIT *ad, *bd;
  1090.    struct descrip td;
  1091.    char tdigits[INTBIGBLK];
  1092.  
  1093.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  1094.       alen = LEN(BigNum(da));
  1095.       blen = LEN(BigNum(db));
  1096.       xlen = alen > blen ? alen : blen;
  1097.       if (blkreq(3 * BigNeed(xlen)) == Error)
  1098.          return Error;
  1099.       a = BigNum(da);
  1100.       b = BigNum(db);
  1101.       x = alcbignum(xlen);
  1102.  
  1103.       if (alen == xlen && !a->sign)
  1104.          ad = DIG(a,0);
  1105.       else {
  1106.          tad = alcbignum(xlen);
  1107.          ad = DIG(tad,0);
  1108.          bdzero(ad, xlen - alen);
  1109.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  1110.          if (a->sign)
  1111.         compl1(ad, ad, xlen);
  1112.          }
  1113.  
  1114.       if (blen == xlen && !b->sign)
  1115.          bd = DIG(b,0);
  1116.       else {
  1117.          tbd = alcbignum(xlen);
  1118.          bd = DIG(tbd,0);
  1119.          bdzero(bd, xlen - blen);
  1120.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  1121.          if (b->sign)
  1122.         compl1(bd, bd, xlen);
  1123.          }
  1124.  
  1125.       for (i = 0; i < xlen; i++)
  1126.          *DIG(x,i) = ad[i] ^ bd[i];
  1127.  
  1128.       if (a->sign ^ b->sign) {
  1129.          x->sign = 1;
  1130.          compl1(DIG(x,0), DIG(x,0), xlen);
  1131.          }
  1132.       }
  1133.    else if (Type(*da) == T_Bignum) {   /* ixor(bignum,integer) */
  1134.       itobig(IntVal(*db), (struct b_bignum *)tdigits, &td);
  1135.       alen = LEN(BigNum(da));
  1136.       blen = LEN(BigNum(&td));
  1137.       xlen = alen > blen ? alen : blen;
  1138.       if (blkreq(3 * BigNeed(alen)) == Error)
  1139.          return Error;
  1140.       a = BigNum(da);
  1141.       b = BigNum(&td);
  1142.       x = alcbignum(alen);
  1143.  
  1144.       if (alen == xlen && !a->sign)
  1145.          ad = DIG(a,0);
  1146.       else {
  1147.          tad = alcbignum(xlen);
  1148.          ad = DIG(tad,0);
  1149.          bdzero(ad, xlen - alen);
  1150.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  1151.          if (a->sign)
  1152.         compl1(ad, ad, xlen);
  1153.          }
  1154.  
  1155.       if (blen == xlen && !b->sign)
  1156.          bd = DIG(b,0);
  1157.       else {
  1158.          tbd = alcbignum(xlen);
  1159.          bd = DIG(tbd,0);
  1160.          bdzero(bd, xlen - blen);
  1161.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  1162.          if (b->sign)
  1163.         compl1(bd, bd, xlen);
  1164.          }
  1165.  
  1166.       for (i = 0; i < xlen; i++)
  1167.          *DIG(x,i) = ad[i] ^ bd[i];
  1168.  
  1169.       if (a->sign ^ b->sign) {
  1170.          x->sign = 1;
  1171.          compl1(DIG(x,0), DIG(x,0), xlen);
  1172.          }
  1173.       }
  1174.    else if (Type(*db) == T_Bignum) {   /* ixor(integer,bignum) */
  1175.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  1176.       alen = LEN(BigNum(&td));
  1177.       blen = LEN(BigNum(db));
  1178.       xlen = alen > blen ? alen : blen;
  1179.       if (blkreq(3 * BigNeed(blen)) == Error)
  1180.          return Error;
  1181.       a = BigNum(&td);
  1182.       b = BigNum(db);
  1183.       x = alcbignum(blen);
  1184.  
  1185.       if (alen == xlen && !a->sign)
  1186.          ad = DIG(a,0);
  1187.       else {
  1188.          tad = alcbignum(xlen);
  1189.          ad = DIG(tad,0);
  1190.          bdzero(ad, xlen - alen);
  1191.          bdcopy(DIG(a,0), &ad[xlen-alen], alen);
  1192.          if (a->sign)
  1193.         compl1(ad, ad, xlen);
  1194.          }
  1195.  
  1196.       if (blen == xlen && !b->sign)
  1197.          bd = DIG(b,0);
  1198.       else {
  1199.          tbd = alcbignum(xlen);
  1200.          bd = DIG(tbd,0);
  1201.          bdzero(bd, xlen - blen);
  1202.          bdcopy(DIG(b,0), &bd[xlen-blen], blen);
  1203.          if (b->sign)
  1204.         compl1(bd, bd, xlen);
  1205.          }
  1206.  
  1207.       for (i = 0; i < xlen; i++)
  1208.          *DIG(x,i) = ad[i] ^ bd[i];
  1209.  
  1210.       if (a->sign ^ b->sign) {
  1211.          x->sign = 1;
  1212.          compl1(DIG(x,0), DIG(x,0), xlen);
  1213.          }
  1214.       }
  1215.    /* not called for ixor(integer,integer) */
  1216.  
  1217.    return mkdesc(x, dx);
  1218. }
  1219.  
  1220. /*
  1221.  *  bigshift(da, db) -> dx
  1222.  */
  1223.  
  1224. int bigshift(da, db, dx)
  1225. dptr da, db, dx;
  1226. {
  1227.    struct b_bignum *a, *x;
  1228.    word alen;
  1229.    word r = IntVal(*db) % NB;
  1230.    word q = (r >= 0 ? IntVal(*db) : (IntVal(*db) - (r += NB))) / NB;
  1231.    word xlen;
  1232.    DIGIT *ad;
  1233.    struct b_bignum *tad;
  1234.    struct descrip td;
  1235.    char tdigits[INTBIGBLK];
  1236.  
  1237.    if (Type(*da) == T_Integer) {
  1238.       itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
  1239.       da = &td;
  1240.       }
  1241.  
  1242.    alen = LEN(BigNum(da));
  1243.    xlen = alen + q + 1;
  1244.    if (xlen <= 0) {
  1245.       MakeInt(0, dx);
  1246.       return Success;
  1247.       }
  1248.    else {
  1249.       if (blkreq(BigNeed(xlen) + BigNeed(alen)) == Error)
  1250.          return Error;
  1251.       a = BigNum(da);
  1252.       x = alcbignum(xlen);
  1253.  
  1254.       if (a->sign && IntVal(*db) > (word)0) {
  1255.          tad = alcbignum(alen);
  1256.          ad = DIG(tad,0);
  1257.          bdcopy(DIG(a,0), ad, alen);
  1258.          compl1(ad, ad, alen);
  1259.          }
  1260.       else
  1261.          ad = DIG(a,0);
  1262.  
  1263.       if (q >= 0) {
  1264.          *DIG(x,0) = shifti1(ad, r, (DIGIT)0, DIG(x,1), alen);
  1265.          bdzero(DIG(x,alen+1), q);
  1266.          }
  1267.       else
  1268.          *DIG(x,0) = shifti1(ad, r, ad[alen+q] >> (NB-r), DIG(x,1), alen+q);
  1269.  
  1270.       if (a->sign && IntVal(*db) > (word)0) {
  1271.          x->sign = 1;
  1272.          *DIG(x,0) |= B - (1 << r);
  1273.          compl1(DIG(x,0), DIG(x,0), xlen);
  1274.          }
  1275.       return mkdesc(x, dx);
  1276.       }
  1277.    }
  1278.  
  1279. /*
  1280.  *  negative if da < db
  1281.  *  zero if da == db
  1282.  *  positive if da > db
  1283.  */
  1284.  
  1285. word bigcmp(da, db)
  1286. dptr da, db;
  1287. {
  1288.    struct b_bignum *a = BigNum(da);
  1289.    struct b_bignum *b = BigNum(db);
  1290.    word alen, blen; 
  1291.  
  1292.    if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
  1293.       if (a->sign != b->sign)
  1294.          return (b->sign - a->sign);
  1295.       alen = LEN(a);
  1296.       blen = LEN(b);
  1297.       if (alen != blen)
  1298.          return (a->sign ? blen - alen : alen - blen);
  1299.  
  1300.       if (a->sign)
  1301.          return cmp1(DIG(b,0), DIG(a,0), alen);
  1302.       else
  1303.          return cmp1(DIG(a,0), DIG(b,0), alen);
  1304.       }
  1305.    else if (Type(*da) == T_Bignum)    /* cmp(bignum, integer) */
  1306.       return bigcmpi(da, IntVal(*db));
  1307.    else                               /* cmp(integer, bignum) */
  1308.       return -bigcmpi(db, IntVal(*da));
  1309. }
  1310.  
  1311. /*
  1312.  *  ?da -> dx
  1313.  */  
  1314.  
  1315. int bigrand(da, dx)
  1316. dptr da, dx;
  1317. {
  1318.    struct b_bignum *a;
  1319.    word alen = LEN(BigNum(da));
  1320.    struct b_bignum *x;
  1321.    struct b_bignum *td;
  1322.    DIGIT *d;
  1323.    word i;
  1324.    double rval;
  1325.  
  1326.    if (blkreq(4 * BigNeed(alen + 1) + 4) == Error)
  1327.       return Error;
  1328.    x = alcbignum(alen);
  1329.    td = alcbignum(alen + 1);
  1330.    d = DIG(td,0);
  1331.    a = BigNum(da);
  1332.  
  1333.    for (i = alen; i >= 0; i--) {
  1334.       rval = RandVal;
  1335.       d[i] = rval * B;
  1336.       }
  1337.     
  1338.    div1(d, DIG(a,0), NULL, DIG(x,0), (word)1, alen);
  1339.    addi1(DIG(x,0), (word)1, DIG(x,0), alen);
  1340.    return mkdesc(x, dx);
  1341. }
  1342.  
  1343. /*
  1344.  *  da + i -> dx
  1345.  */
  1346.  
  1347. static int bigaddi(da, i, dx)
  1348. dptr da, dx;
  1349. word i;
  1350. {
  1351.    struct b_bignum *a, *x; 
  1352.    word alen; 
  1353.  
  1354.    if (i < 0)
  1355.       return bigsubi(da, -i, dx);
  1356.    else if (i != (DIGIT)i) {
  1357.       struct descrip td;
  1358.       char tdigits[INTBIGBLK];
  1359.  
  1360.       itobig(i, (struct b_bignum *)tdigits, &td);
  1361.       return bigadd(da, &td, dx);
  1362.       }
  1363.    else {
  1364.       alen = LEN(BigNum(da));
  1365.       if (blkreq(BigNeed(alen + 1)) == Error)
  1366.          return Error;
  1367.       a = BigNum(da);
  1368.       if (a->sign) {
  1369.      x = alcbignum(alen);
  1370.          subi1(DIG(a,0), i, DIG(x,0), alen);
  1371.          }
  1372.       else {
  1373.          x = alcbignum(alen + 1);
  1374.          *DIG(x,0) = addi1(DIG(a,0), i, DIG(x,1), alen);
  1375.          }
  1376.       x->sign = a->sign;
  1377.       return mkdesc(x, dx);
  1378.       }
  1379. }
  1380.  
  1381. /*
  1382.  *  da - i -> dx
  1383.  */
  1384.  
  1385. static int bigsubi(da, i, dx)
  1386. dptr da, dx;
  1387. word i;
  1388. {
  1389.    struct b_bignum *a, *x; 
  1390.    word alen;
  1391.  
  1392.    if (i < 0)
  1393.       return bigaddi(da, -i, dx);
  1394.    else if (i != (DIGIT)i) {
  1395.       struct descrip td;
  1396.       char tdigits[INTBIGBLK];
  1397.  
  1398.       itobig(i, (struct b_bignum *)tdigits, &td);
  1399.       return bigsub(da, &td, dx);
  1400.       }
  1401.    else {
  1402.       alen = LEN(BigNum(da));
  1403.       if (blkreq(BigNeed(alen + 1)) == Error)
  1404.      return Error;
  1405.       a = BigNum(da);
  1406.       if (a->sign) {
  1407.          x = alcbignum(alen + 1);
  1408.          *DIG(x,0) = addi1(DIG(a,0), i, DIG(x,1), alen);
  1409.          }
  1410.       else {
  1411.          x = alcbignum(alen);
  1412.          subi1(DIG(a,0), i, DIG(x,0), alen);
  1413.          }
  1414.       x->sign = a->sign;
  1415.       return mkdesc(x, dx);
  1416.       }
  1417. }
  1418.  
  1419. /*
  1420.  *  da * i -> dx
  1421.  */
  1422.  
  1423. static int bigmuli(da, i, dx)
  1424. dptr da, dx;
  1425. word i;
  1426. {
  1427.    struct b_bignum *a, *x; 
  1428.    word alen;
  1429.  
  1430.    if (i <= -B || i >= B) {
  1431.       struct descrip td;
  1432.       char tdigits[INTBIGBLK];
  1433.  
  1434.       itobig(i, (struct b_bignum *)tdigits, &td);
  1435.       return bigmul(da, &td, dx);
  1436.       }
  1437.    else {
  1438.       alen = LEN(BigNum(da));
  1439.       if (blkreq(BigNeed(alen + 1)) == Error)
  1440.      return Error;
  1441.       a = BigNum(da);
  1442.       x = alcbignum(alen + 1);
  1443.       if (i >= 0)
  1444.          x->sign = a->sign;
  1445.       else {
  1446.          x->sign = 1 ^ a->sign;
  1447.          i = -i;
  1448.          }
  1449.       *DIG(x,0) = muli1(DIG(a,0), i, 0, DIG(x,1), alen);
  1450.       return mkdesc(x, dx);
  1451.       }
  1452. }
  1453.  
  1454. /*
  1455.  *  da / i -> dx
  1456.  */
  1457.  
  1458. static int bigdivi(da, i, dx)
  1459. dptr da, dx;
  1460. word i;
  1461. {
  1462.    struct b_bignum *a, *x; 
  1463.    word alen;
  1464.  
  1465.    if (i <= -B || i >= B) {
  1466.       struct descrip td;
  1467.       char tdigits[INTBIGBLK];
  1468.  
  1469.       itobig(i, (struct b_bignum *)tdigits, &td);
  1470.       return bigdiv(da, &td, dx);
  1471.       }
  1472.    else {
  1473.       alen = LEN(BigNum(da));
  1474.       if (blkreq(BigNeed(alen)) == Error)
  1475.      return Error;
  1476.       a = BigNum(da);
  1477.       x = alcbignum(alen);
  1478.       if (i >= 0)
  1479.          x->sign = a->sign;
  1480.       else {
  1481.          x->sign = 1 ^ a->sign;
  1482.          i = -i;
  1483.          }
  1484.       divi1(DIG(a,0), i, DIG(x,0), alen);
  1485.       return mkdesc(x, dx);
  1486.       }
  1487. }
  1488.  
  1489. /*
  1490.  *  da % i -> dx
  1491.  */
  1492.  
  1493. static int bigmodi(da, i, dx)
  1494. dptr da, dx;
  1495. word i;
  1496. {
  1497.    struct b_bignum *a, *temp;
  1498.    word alen;
  1499.    word x;
  1500.  
  1501.    if (i <= -B || i >= B) {
  1502.       struct descrip td;
  1503.       char tdigits[INTBIGBLK];
  1504.  
  1505.       itobig(i, (struct b_bignum *)tdigits, &td);
  1506.       return bigmod(da, &td, dx);
  1507.       }
  1508.    else {
  1509.       alen = LEN(BigNum(da));
  1510.       if (blkreq(BigNeed(alen)) == Error)
  1511.      return Error;
  1512.       a = BigNum(da);
  1513.       temp = alcbignum(alen);
  1514.       x = divi1(DIG(a,0), Abs(i), DIG(temp,0), alen);
  1515.       if (a->sign)
  1516.      x = -x;
  1517.       MakeInt(x, dx);
  1518.       return Success;
  1519.       }
  1520. }
  1521.  
  1522. /*
  1523.  *  da ^ i -> dx
  1524.  */
  1525.  
  1526. static int bigpowi(da, i, dx)
  1527. dptr da, dx;
  1528. word i;
  1529. {
  1530.    int n = WordBits;
  1531.  
  1532.    if (i > 0) {
  1533.       /* scan bits left to right.  skip leading 1. */
  1534.       while (--n >= 0)
  1535.          if (i & ((word)1 << n))
  1536.         break;
  1537.       /* then, for each zero, square the partial result;
  1538.          for each one, square it and multiply it by a */
  1539.       *dx = *da;
  1540.       while (--n >= 0) {
  1541.          if (bigmul(dx, dx, dx) == Error)
  1542.         return Error;
  1543.          if (i & ((word)1 << n))
  1544.             if (bigmul(dx, da, dx) == Error)
  1545.            return Error;
  1546.          }
  1547.       }
  1548.    else {
  1549.       MakeInt(0, dx);
  1550.       }
  1551.    return Success;
  1552. }
  1553.  
  1554. /*
  1555.  *  a ^ i -> dx
  1556.  */
  1557.  
  1558. static int bigpowii(a, i, dx)
  1559. word a, i;
  1560. dptr dx;
  1561. {
  1562.    word x, y;
  1563.    int n = WordBits;
  1564.    int isbig = 0;
  1565.  
  1566.    if (a == 0 || i <= 0) {              /* special cases */
  1567.       if (a == 0 && i <= 0)             /* 0 ^ negative -> error */
  1568.          RetError(-204, nulldesc);
  1569.       if (a == -1) {                    /* -1 ^ [odd,even] -> [-1,+1] */
  1570.          if (!(i & 1))
  1571.         a = 1;
  1572.          }
  1573.       else if (a != 1) {                /* 1 ^ any -> 1 */
  1574.          a = 0;
  1575.          }                   /* others ^ negative -> 0 */
  1576.       MakeInt(a, dx);
  1577.       }
  1578.    else {
  1579.       struct descrip td;
  1580.       char tdigits[INTBIGBLK];
  1581.  
  1582.       /* scan bits left to right.  skip leading 1. */
  1583.       while (--n >= 0)
  1584.          if (i & ((word)1 << n))
  1585.         break;
  1586.       /* then, for each zero, square the partial result;
  1587.          for each one, square it and multiply it by a */
  1588.       x = a;
  1589.       while (--n >= 0) {
  1590.          if (isbig) {
  1591.             if (bigmul(dx, dx, dx) == Error)
  1592.            return Error;
  1593.         }
  1594.          else {
  1595.             y = mul(x, x);
  1596.             if (!over_flow)
  1597.                x = y;
  1598.             else {
  1599.                itobig(x, (struct b_bignum *)tdigits, &td);
  1600.                if (bigmul(&td, &td, dx) == Error)
  1601.                  return Error;
  1602.                isbig = (Type(*dx) == T_Bignum);
  1603.                } 
  1604.             }
  1605.          if (i & ((word)1 << n)) {
  1606.             if (isbig) {
  1607.                if (bigmuli(dx, a, dx) == Error)
  1608.           return Error;
  1609.            }
  1610.             else {
  1611.                y = mul(x, a);
  1612.                if (!over_flow)
  1613.                   x = y;
  1614.                else {
  1615.                   itobig(x, (struct b_bignum *)tdigits, &td);
  1616.                   if (bigmuli(&td, a, dx) == Error)
  1617.              return Error;
  1618.                   isbig = (Type(*dx) == T_Bignum);
  1619.                   }
  1620.                }
  1621.             }
  1622.          }
  1623.       if (!isbig) {
  1624.      MakeInt(x, dx);
  1625.      }
  1626.       }
  1627.    return Success;
  1628. }
  1629.  
  1630. /*
  1631.  *  negative if da < i
  1632.  *  zero if da == i
  1633.  *  positive if da > i
  1634.  */  
  1635.   
  1636. static word bigcmpi(da, i)
  1637. dptr da;
  1638. word i;
  1639. {
  1640.    struct b_bignum *a = BigNum(da);
  1641.    word alen = LEN(a);
  1642.  
  1643.    if (i > -B && i < B) {
  1644.       if (i >= 0)
  1645.          if (a->sign)
  1646.         return -1;
  1647.          else
  1648.         return cmpi1(DIG(a,0), i, alen);
  1649.       else
  1650.          if (a->sign)
  1651.         return -cmpi1(DIG(a,0), -i, alen);
  1652.          else
  1653.         return 1;
  1654.       }
  1655.    else {
  1656.       struct descrip td;
  1657.       char tdigits[INTBIGBLK];
  1658.  
  1659.       itobig(i, (struct b_bignum *)tdigits, &td);
  1660.       return bigcmp(da, &td);
  1661.       }
  1662. }
  1663.  
  1664.  
  1665. /* These are all straight out of Knuth vol. 2, Sec. 4.3.1. */
  1666.  
  1667. /*
  1668.  *  (u,n) + (v,n) -> (w,n)
  1669.  *
  1670.  *  returns carry, 0 or 1
  1671.  */
  1672.  
  1673. static DIGIT add1(u, v, w, n)
  1674. DIGIT *u, *v, *w;
  1675. word n;
  1676. {
  1677.    uword dig, carry; 
  1678.    word i;
  1679.  
  1680.    carry = 0;
  1681.    for (i = n; --i >= 0; ) {
  1682.       dig = (uword)u[i] + v[i] + carry;
  1683.       w[i] = lo(dig);
  1684.       carry = hi(dig);
  1685.       }
  1686.    return carry;
  1687. }
  1688.  
  1689. /*
  1690.  *  (u,n) - (v,n) -> (w,n)
  1691.  *
  1692.  *  returns carry, 0 or -1
  1693.  */
  1694.  
  1695. static word sub1(u, v, w, n)
  1696. DIGIT *u, *v, *w;
  1697. word n;
  1698. {
  1699.    uword dig, carry; 
  1700.    word i;
  1701.  
  1702.    carry = 0;
  1703.    for (i = n; --i >= 0; ) {
  1704.       dig = (uword)u[i] - v[i] + carry;
  1705.       w[i] = lo(dig);
  1706.       carry = signed_hi(dig);
  1707.       }
  1708.    return carry;
  1709. }
  1710.  
  1711. /*
  1712.  *  (u,n) * (v,m) -> (w,m+n)
  1713.  */
  1714.  
  1715. static novalue mul1(u, v, w, n, m)
  1716. DIGIT *u, *v, *w;
  1717. word n, m;
  1718. {
  1719.    word i, j;
  1720.    uword dig, carry;
  1721.  
  1722.    bdzero(&w[m], n);
  1723.  
  1724.    for (j = m; --j >= 0; ) {
  1725.       carry = 0;
  1726.       for (i = n; --i >= 0; ) {
  1727.          dig = (uword)u[i] * v[j] + w[i+j+1] + carry;
  1728.          w[i+j+1] = lo(dig);
  1729.          carry = hi(dig);
  1730.          }
  1731.       w[j] = carry;
  1732.       }
  1733. }
  1734.  
  1735. /*
  1736.  *  (a,m+n) / (b,n) -> (q,m+1) (r,n)
  1737.  *
  1738.  *  if q or r is NULL, the quotient or remainder is discarded
  1739.  */
  1740.  
  1741. static novalue div1(a, b, q, r, m, n)
  1742. DIGIT *a, *b, *q, *r;
  1743. word m, n;
  1744. {
  1745.    uword qhat, rhat;
  1746.    uword dig, carry;
  1747.    struct b_bignum *tu, *tv;
  1748.    DIGIT *u, *v;
  1749.    word d;
  1750.    word i, j;
  1751.  
  1752.    /* blkreq's done in calling routines */
  1753.  
  1754.    tu = alcbignum(m + n + 1);
  1755.    tv = alcbignum(n);
  1756.    u = DIG(tu,0);
  1757.    v = DIG(tv,0);
  1758.  
  1759.    /* D1 */
  1760.    for (d = 0; d < NB; d++)
  1761.       if (b[0] & (1 << (NB - 1 - d)))
  1762.          break;
  1763.  
  1764.    u[0] = shifti1(a, d, (DIGIT)0, &u[1], m+n);
  1765.    shifti1(b, d, (DIGIT)0, v, n);
  1766.  
  1767.    /* D2, D7 */
  1768.    for (j = 0; j <= m; j++) {
  1769.       /* D3 */
  1770.       if (u[j] == v[0]) {
  1771.          qhat = B - 1;
  1772.          rhat = (uword)v[0] + u[j+1];
  1773.          }
  1774.       else {
  1775.          uword numerator = dbl(u[j], u[j+1]);
  1776.          qhat = numerator / (uword)v[0];
  1777.          rhat = numerator % (uword)v[0];
  1778.          }
  1779.  
  1780.       while (rhat < B && qhat * v[1] > dbl(rhat, u[j+2])) {
  1781.          qhat -= 1;
  1782.          rhat += v[0];
  1783.          }
  1784.             
  1785.       /* D4 */
  1786.       carry = 0;
  1787.       for (i = n; i > 0; i--) {
  1788.          dig = u[i+j] - v[i-1] * qhat + carry;       /* -BSQ+B .. B-1 */
  1789.          u[i+j] = lo(dig);
  1790.          if ((uword)dig < B)
  1791.             carry = hi(dig);
  1792.          else carry = hi(dig) | -B;
  1793.          }
  1794.       carry = (word)(carry + u[j]) < 0;
  1795.  
  1796.       /* D5 */
  1797.       if (q)
  1798.      q[j] = qhat;
  1799.  
  1800.       /* D6 */
  1801.       if (carry) {
  1802.          if (q)
  1803.         q[j] -= 1;
  1804.          carry = 0;
  1805.          for (i = n; i > 0; i--) {
  1806.             dig = (uword)u[i+j] + v[i-1] + carry;
  1807.             u[i+j] = lo(dig);
  1808.             carry = hi(dig);
  1809.             }
  1810.          }
  1811.       }
  1812.  
  1813.    if (r) {
  1814.       if (d == 0)
  1815.          shifti1(&u[m+1], (word)d, (DIGIT)0, r, n);
  1816.       else
  1817.          r[0] = shifti1(&u[m+1], (word)(NB - d), u[m+n]>>d, &r[1], n - 1);
  1818.       }
  1819. }
  1820.  
  1821. /*
  1822.  *  - (u,n) -> (w,n)
  1823.  *
  1824.  */
  1825.  
  1826. static novalue compl1(u, w, n)
  1827. DIGIT *u, *w;
  1828. word n;
  1829. {
  1830.    uword dig, carry = 0;
  1831.    word i;
  1832.  
  1833.    for (i = n; --i >= 0; ) {
  1834.       dig = carry - u[i];
  1835.       w[i] = lo(dig);
  1836.       carry = signed_hi(dig);
  1837.       }
  1838. }
  1839.  
  1840. /*
  1841.  *  (u,n) : (v,n)
  1842.  */
  1843.  
  1844. static word cmp1(u, v, n)
  1845. DIGIT *u, *v;
  1846. word n;
  1847. {
  1848.    word i;
  1849.  
  1850.    for (i = 0; i < n; i++)
  1851.       if (u[i] != v[i])
  1852.          return u[i] > v[i] ? 1 : -1;
  1853.    return 0;
  1854. }
  1855.  
  1856. /*
  1857.  *  (u,n) + k -> (w,n)
  1858.  *
  1859.  *  k in 0 .. B-1
  1860.  *  returns carry, 0 or 1
  1861.  */
  1862.  
  1863. static DIGIT addi1(u, k, w, n)
  1864. DIGIT *u, *w;
  1865. word k;
  1866. word n;
  1867. {
  1868.    uword dig, carry;
  1869.    word i;
  1870.     
  1871.    carry = k;
  1872.    for (i = n; --i >= 0; ) {
  1873.       dig = (uword)u[i] + carry;
  1874.       w[i] = lo(dig);
  1875.       carry = hi(dig);
  1876.       }
  1877.    return carry;
  1878. }
  1879.  
  1880. /*
  1881.  *  (u,n) - k -> (w,n)
  1882.  *
  1883.  *  k in 0 .. B-1
  1884.  *  u must be greater than k
  1885.  */
  1886.  
  1887. static novalue subi1(u, k, w, n)
  1888. DIGIT *u, *w;
  1889. word k;
  1890. word n;
  1891. {
  1892.    uword dig, carry;
  1893.    word i;
  1894.     
  1895.    carry = -k;
  1896.    for (i = n; --i >= 0; ) {
  1897.       dig = (uword)u[i] + carry;
  1898.       w[i] = lo(dig);
  1899.       carry = signed_hi(dig);
  1900.       }
  1901. }
  1902.  
  1903. /*
  1904.  *  (u,n) * k + c -> (w,n)
  1905.  *
  1906.  *  k in 0 .. B-1
  1907.  *  returns carry, 0 .. B-1
  1908.  */
  1909.  
  1910. static DIGIT muli1(u, k, c, w, n)
  1911. DIGIT *u, *w;
  1912. word k;
  1913. int c;
  1914. word n;
  1915. {
  1916.    uword dig, carry;
  1917.    word i;
  1918.  
  1919.    carry = c;
  1920.    for (i = n; --i >= 0; ) {
  1921.       dig = (uword)k * u[i] + carry;
  1922.       w[i] = lo(dig);
  1923.       carry = hi(dig);
  1924.       }
  1925.    return carry;
  1926. }
  1927.  
  1928. /*
  1929.  *  (u,n) / k -> (w,n)
  1930.  *
  1931.  *  k in 0 .. B-1
  1932.  *  returns remainder, 0 .. B-1
  1933.  */
  1934.  
  1935. static DIGIT divi1(u, k, w, n)
  1936. DIGIT *u, *w;
  1937. word k;
  1938. word n;
  1939. {
  1940.    uword dig, remain;
  1941.    word i;
  1942.  
  1943.    remain = 0;
  1944.    for (i = 0; i < n; i++) {
  1945.       dig = dbl(remain, u[i]);
  1946.       w[i] = dig / k;
  1947.       remain = dig % k;
  1948.       }
  1949.    return remain;
  1950. }
  1951.  
  1952. /*
  1953.  *  ((u,n) << k) + c -> (w,n)
  1954.  *
  1955.  *  k in 0 .. NB-1
  1956.  *  c in 0 .. B-1 
  1957.  *  returns carry, 0 .. B-1
  1958.  */
  1959.  
  1960. static DIGIT shifti1(u, k, c, w, n)
  1961. DIGIT *u, c, *w;
  1962. word k;
  1963. word n;
  1964. {
  1965.    uword dig;
  1966.    word i;
  1967.  
  1968.    if (k == 0) {
  1969.       bdcopy(u, w, n);
  1970.       return 0;
  1971.       }
  1972.     
  1973.    for (i = n; --i >= 0; ) {
  1974.       dig = ((uword)u[i] << k) + c;
  1975.       w[i] = lo(dig);
  1976.       c = hi(dig);
  1977.       }
  1978.    return c;
  1979. }
  1980.  
  1981. /*
  1982.  *  (u,n) : k
  1983.  *
  1984.  *  k in 0 .. B-1
  1985.  */
  1986.  
  1987. static word cmpi1(u, k, n)
  1988. DIGIT *u;
  1989. word k;
  1990. word n;
  1991. {
  1992.    word i;
  1993.  
  1994.    for (i = 0; i < n-1; i++)
  1995.       if (u[i])
  1996.      return 1;
  1997.    if (u[n - 1] == (DIGIT)k)
  1998.       return 0;
  1999.    return u[n - 1] > (DIGIT)k ? 1 : -1;
  2000. }
  2001.  
  2002. static novalue bdzero(dest, l)
  2003. DIGIT *dest;
  2004. word l;
  2005. {
  2006.    word i;
  2007.  
  2008.    for (i = 0; i < l; i++)
  2009.       dest[i] = 0;
  2010. }
  2011.  
  2012. static novalue bdcopy(src, dest, l)
  2013. DIGIT *src, *dest;
  2014. word l;
  2015. {
  2016.    word i;
  2017.  
  2018.    for (i = 0; i < l; i++)
  2019.       dest[i] = src[i];
  2020. }
  2021. #else                    /* LargeInts */
  2022. static char x;            /* prevent empty module */
  2023. #endif                    /* LargeInts */
  2024.